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^^Jie Density Matrix Renormalization Group (DMRG) has become a powerful numerical method that 
jl^an be applied to low-dimensional strongly correlated fermionic and bosonic systems. It allows for 
■^J^arvery precise calculation of static, dynamic and thermodynamic properties. Its field of applicabil- 
ity has now extended beyond Condensed Matter, and is successfully used in Quantum Chemistry, 
'+j5tatistical Mechanics, Quantum Information Theory, Nuclear and High Energy Physics as well. In 
Tjhis article, we briefly review the main aspects of the method and present some of the most rele- 
^-Hvant applications so as to give an overview on the scope and possibilities of DMRG. We focus on 
I the most important extensions of the method such as the calculation of dynamical properties, the 
'"^^application to classical systems, finite temperature simulations, phonons and disorder, field theory, 
(-Htime-dependent properties and the ab initio calculation of electronic states in molecules. The recent 
/-sjquantum information interpretation, the development of highly accurate time-dependent algorithms 
i^and the possibility of using the DMRG as the impurity-solver of the Dynamical Mean Field Method 
I '•^CDMFT') give new insights into its present and potential uses. We review the numerous very recent 
applications of these techniques where the DMRG has shown to be one of the most reliable and 
^^rsa,tile methods in modern computational physics. 
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1 Introduction 

The Density Matrix Renormalization Group (DMRG) was developed by S. 
White in 1992 and since then it has proved to be a very powerful method 
for low dimensional interacting systems. Its remarkable accuracy can be seen 
for example in the spin-1 Heisenberg chain: for a system of hundreds of sites 
a precision of 10"^'' for the ground state energy can be achieved. Since then it 
has been applied to a great variety of systems and problems including, among 
others, spin chains and ladders, fermionic and bosonic systems, disordered 
models, impurities, molecules, nanoscopic systems and 2D electrons in high 
magnetic fields. It has also been improved substantially in several directions 
like two dimensional (2D) classical systems, stochastic models, inclusion of 
phonons, quantum chemistry, field theory, finite temperature and the calcula- 
tion of dynamical and time-dependent properties. Some calculations have also 
been performed in 2D quantum systems. Most of these topics are treated in 
detail and in a pedagogical way in the book where the reader can find an 
extensive overview on the foundations of DMRG and also in a comprehensive 
review Recent new developments of the DMRG like the implementation 
of very accurate methods for time-dependent problems, the quantum informa- 
tion and matrix-product perspective and the possibility to combine it with the 
Dynamical Mean Field Theory (DMFT) , triggered a great activity using these 
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techniques, as seen by the numerous recent papers pubhshed in the last year. 
In this article we will mainly focus on these new developments and applica- 
tions and hope to give new insights into its potential uses. In order to achieve 
this and on behalf of coherence we will also briefly describe the method and 
its former extensions. 

1.1 Motivation 

When considering finite systems, the exponential growth of degrees of freedom 
to be considered imposes an important limitation in numerical calculations. 
Several methods have been introduced in order to reduce the size of the Hilbert 
space to be able to reach larger systems, such as Monte Carlo, renormalization 
group (RG) and DMRG. Each method considers a particular criterion for 
keeping the relevant information. 

The DMRG was originally developed to overcome the problems that arise 
in interacting systems in ID when standard RG procedures were applied. For 
example, consider a block B (a block is a collection of sites) where the Hamil- 
tonian Hb and end-operators are defined. These traditional methods consist 
in putting together two or more blocks (e.g. B-B', which we will call the su- 
perblock), connected using end-operators, in a basis that is a direct product of 
the basis of each block, forming Hbb'- This Hamiltonian is then diagonalized, 
the superblock is replaced by a new effective block Bnew formed by a certain 
number m of lowest-lying eigenstates of Hbb' and the iteration is continued 
(see Ref. Although it has been used successfully in certain cases, this pro- 
cedure, or similar versions of it, has been applied to several interacting systems 
with poor performance. For example, it has been applied to the ID Hubbard 
model keeping m ~ 1000 states and for 16 sites, an error of 5-10% was ob- 
tained j^. Other results ^ were also discouraging. A better performance was 
obtained by adding a single site at a time rather than doubling the block 
size. However, there is one case where a similar version of this method ap- 
plies very well: the impurity problem represented by the Kondo and Anderson 
models. Wilson |H] mapped the one-impurity problem onto a one-dimensional 
lattice with exponentially decreasing hoppings. The difference with the method 
explained above is that in this case, one site (equivalent to an "onion shell") 
is added at each step and, due to the exponential decrease of the hopping, 
very accurate results can be obtained. A very recent work on renormalization 
group transformations on quantum states based on matrix product states was 
performed in Ref. |9|. 

Returning to the problem of putting several blocks together, the main source 
of error comes from the election of eigenstates of Hbb' as representative states 
of a superblock. Since Hbb' has no connection to the rest of the lattice, its 
eigenstates may have unwanted features (like nodes) at the ends of the block 
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and this can't be improved by increasing the number of states kept. Based on 
this consideration, Noack and White JTI3|^ tried including different boundary 
conditions and boundary strengths. This turned out to work well for single 
particle and Anderson localization problems but, however, it did not improve 
the results significantly for interacting systems. These considerations led to 
the idea of taking a larger superblock that includes the blocks BB\ diagonal- 
ize the Hamiltonian in this large superblock and then somehow project the 
most favorable states onto BB' . Then BB' is replaced by Bnew In this way, 
awkward features in the boundary would vanish and a better representation 
of the states in the infinite system would be achieved. White appealed to 
Feynman's formulation of the density matrix as the best description of a part 
of a quantum mechanical system and, thus, as the optimal way of projecting 
the most relevant states onto a subsystem. In Ref. jllj . Gaite shows that, by 
using an angular quantization construction of the density matrix, the DMRG 
is an algorithm that keeps states with a higher weight near the boundaries 
in a systematic way. This is a consequence of the fact that, near the bound- 
aries, there is a concentration of quantum states, in a similar manner as in the 
physics of black holes. 

In the following Section we will describe the standard method; in Sect. 3 we 
will mention some of the most important applications, describing in Sect. 4 the 
treatment of bosonic degrees of freedom. In Sec. 5 we review the most relevant 
extensions of the method like the momentum and energy levels representation 
and the application to Quantum Chemistry. The improvement to treat systems 
in dimensions higher than 1 is mentioned in Sec. 6, while Sect. 7 will be devoted 
to the latest and very promising developments to calculate non-equilibrium 
and time-dependent properties. In Sect. 8 we deal with the calculation of 
dynamical properties at zero temperature and in Sect. 9 we briefly describe 
how DMRG can be used to calculate physical properties of classical systems. 
In connection to this, finite temperature studies were possible within DMRG 
and this is explained in Section 10. In the following Sect. 11 we describe how 
the Dynamical Mean Field Theory (DMFT) can profit from the DMRG as an 
impurity-solver to obtain spectral properties of correlated systems, opening 
new posibilities for the calculation of more complicated and realistic systems. 
Finally, we summarize the achievements of the method and look upon new 
potential applications and developments. 



2 The DMRG Method 

The DMRG allows for a systematic truncation of the Hilbert space by keeping 
the most probable states describing a wave function {e.g. the ground state) 
instead of the lowest energy states usually kept in previous real space renor- 
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System Environment 




Figure 1. A scheme of the super block (universe) configuration for the DMRG algorithm |4|. 

malization techniques. 

The basic idea consists in starting with a smah system {e.g with sites) 
and then gradually increase its size (to A^+4,...) until the desired length 

is reached. Let us call the collection of N sites the universe and divide it into 
two parts: the system and the environment (see Fig. 1). The Hamiltonian 
is constructed in the universe and its ground state |'0o) is obtained. This 
is considered as the state of the universe and called the target state. It has 
components on the system and the environment. We want to obtain the most 
relevant states of the system, i.e., the states of the system that have largest 
weight in |'i/'o)- To obtain this, the environment is considered as a statistical 
bath and the density matrix jT2] is used to obtain the desired information on 
the system. So instead of keeping eigenstates of the Hamiltonian in the block 
(system), we keep eigenstates of the density matrix. We will be more explicit 
below. 

A very easy and pedagogical way of understanding the basic functioning 
of DMRG is applying it to the calculation of simple quantum problems like 
one particle in a tight-binding chain |13[ ll4j. In these examples, a discretized 
version of the Schrodinger equation is used and extremely accurate results are 
obtained for the quantum harmonic oscillator, the anharmonic oscillator and 
the double-well potential. 

Let's define block [B] as a finite chain with I sites having an associated 
Hilbert space with M total amount of states (or m states after a reduction 
is performed) where operators are defined (in particular the Hamiltonian in 
this finite chain, Hb, and the operators at the ends of the block, useful for 
linking it to other chains or added sites). Except for the first iteration, the 
basis in this block isn't explicitly known due to previous basis rotations and 
reductions. The operators in this basis are matrices and the basis states are 
characterized by quantum numbers (like S^, charge or number of particles, 
etc). We also define an added block or site as [a] having n states. A general 
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iteration of the method is described below: 

i) Define the Hamiltonian Hbb' for the superblock (the universe) formed 
by putting together two blocks [B] and [B'] and two added sites [a] and [a'] 
in this way: [B a a' B' ], where the primes indicate additional blocks. The 
primed blocks usually have the same structure as the non-primed ones, but 
this can vary (see the finite-system algorithm below). The block [B'] has M' 
or m' states for the full and reduced spaces, respectively. In general, blocks 
[B] and [B'] come from the previous iteration. The total Hilbert space of this 
superblock is the direct product of the individual spaces corresponding to each 
block and the added sites. In practice a quantum number of the superblock 
can be fixed (in a spin chain for example one can look at the total = 
subspace), so the total number of states in the superblock is much smaller 
than (mn)^. In some cases, as the quantum number of the superblock consists 
of the sum of the quantum numbers of the individual blocks, each block must 
contain several subspaces (several values of 5"^ for example). Here periodic 
boundary conditions can be attached to the ends and a different block layout 
should be considered (e.g. [B a B' a' ]) to avoid connecting blocks [B] and 
[B'] which takes longer to converge. The boundary conditions are between [a'] 
and [B]. For closed chains the performance is poorer than for open boundary 
conditions [HITK] (see Fig. 2) 

ii) Diagonalize the Hamiltonian Hbb' to obtain the ground state |^o) (target 
state) using Lanczos '16' or Davidson HT^ algorithms. Other states could also 
be kept, such as the first excited ones: they are all called target states. A faster 
convergence of Lanczos or Davidson algorithm is achieved by choosing a good 
trial vector [THllTn] . 

iii) Construct the density matrix: 



on block [B a], where V'o.ij = jiV'o)) the states \i) belonging to the Hilbert 
space of the block [B a] and the states |j) to the block [B' a']. The den- 
sity matrix considers the part [B a] as a system and [B' a'], as a statistical 
bath. The eigenstates of p with the highest eigenvalues correspond to the 
most probable states (or equivalently the states with highest weight) of block 
[B a] in the ground state of the whole superblock. These states are kept up 
to a certain cutoff, keeping a total of m states per block. The density ma- 
trix eigenvalues, uja, sum up to unity and the truncation error, defined as the 
sum of the density matrix eigenvalues corresponding to discarded eigenvectors, 
X]™=m+i = ^ — sum^^iUJa, gives a qualitative indication of the accuracy of 
the calculation {rrimax corresponds to the size of the full space of block [B a]), 
iv) With these m states a rectangular matrix O is formed and it is used to 
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Figure 2. Truncation error for the finite-system DMRG as a function of the number of states kept 
m (reprinted from •13' , with permission). 

change basis and reduce all operators defined in [B a]. This block [B a] is 
then renamed as block [B^gu;] or simply [B] (for example, the Hamiltonian in 
block [B a], HBa, is transformed into Hb as Hb = O'^HsaO). 

v) A new block [a] is added (one site in our case) and the new superblock 
[B a a' B'] is formed as the direct product of the states of all the blocks. 

vi) This iteration continues until the desired length is achieved. At each step 
the length is = 2/ + 2 (if [a] consists of one site). 

When more than one target state is used, i.e more than one state is wished 
to be well described, the density matrix is defined as: 

Pii' = ^Pi^ 4'i,ij4'i,i'j (2) 
I j 

where pi defines the probability of finding the system in the target state {(pi) 
(not necessarily eigenstates of the Hamiltonian) . 

The method described above is usually called the infinite-system algorithm 
since the system size increases at each iteration. There is a way to increase 
precision at each length A^ called the finite-system algorithm. It consists of 
fixing the lattice size and zipping to and fro a couple of times until convergence 
is reached. In this case and for the block configuration [B a a' B' ], A^ = 
/ + 1 + 1 + /' where / and I' are the number of sites in B and B' respectively. In 
this step the density matrix is used to project onto the left / + 1 sites (see Fig. 
121). In order to keep A^ fixed, in the next block configuration, the right block 
B' should be defined in i - 1 sites such that N = {I + 1) + 1 + 1 + {I - 1)'. The 
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Figure 3. One step of the finite-system algorithm, where the total length N = I + I' + 2 is fixed. 
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Figure 4. Difi'erence between the ground-state energy obtained from the finite- system DMRG with 
different number of iterations and number of states kept m, and the exact energy calculated using 
Bethe Ansatz for a 32-site Hubbard model (U = 2 and filling n = 3/4). Reprinted from |13| with 

permission. 

operators in this smaller block should be kept from previous iterations (in some 
cases from the iterations for the system size with N — 2) 2 . With this method, 
a higher precision can be achieved (see Fig. 0]). It is only within this finite- 
system algorithm that the total errors are controlled (and are proportional) 
to the discarded error of the density matrix. 

As an example of a general DMRG building step let's look at the one- 
dimensional Heisenberg model 

-ff = ^ SjSj+i = SfSf^i + -{SfS~^^ + S^Sfj^^). (3) 

i 

Suppose we want to connect two blocks [Si=B a ], defined by states with 
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[i?2=a'B'], defined by states i' ,j' , to form tlie Hamiltonian in tlie superblock, 
wfiicli we will diagonalize to get the energies and eigenstates. Performing an 
external product between the operators we have stored from previous iterations 
we have: 

[HB,B,]ir,i'r = [HB,]ii'Sjj' + [HB,]jfSii, + [S'Jii^ [5^,]^,^ (4) 

+ 2\-^a]ii'[^a']jj' + 2i^a ]ii'[^a']jj' 

If fermionic models are to be considered, more operators should be kept, 
such as, for example, the creation and destruction operators at the borders of 
the blocks to perform the interblock hopping. 

The calculation of static properties like correlation functions is easily done 
by keeping the operators at each step and performing the corresponding basis 
change and reduction, in a similar manner as done with the Hamiltonian in 
each block [3]. The energy and measurements are calculated in the superblock. 

2.1 The density matrix projection 

The density matrix leads naturally to the optimal states in the system as we 
demonstrate below [T^Oj. 
Let 

M 

m=Y.^^Mj) (5) 

be a state of the universe (system + environment), having real coefficients 
for simplicity (here we will assume the same configuration of system and en- 
vironment so M = M' , but this assumption is not necessary). We want to 
obtain a variational wave function \iPq) defined in an optimally reduced space, 
generated by the m "system" vectors \a) = YliLiUai\i) 

m M 

\'<1jo) = ^^o,aM)\j) (6) 

a=l j=l 

such that the modulus of the difference with the true wave functional 



aij aj 
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is minimal w.r.t the a^j. Setting the derivative in these variables to zero, this 
condition leads to 

^Vii^^m = aaj- (8) 
i 

If we define the density operator for a pure state of the "universe" p = |'0o) (V'ol > 
the reduced density matrix of the "system" is: 

Pii' = ^Pij,i'j = ^{j\{i\i'o){i^o\i')\j) = ^i'iji'i'j, (9) 

j 3 3 

Replacing |H1 in Eq. [3 and using the above expression we obtain the condition: 

m 

1 - ^ UaiPii'Uai' = 1 - ^ Wa, (10) 
aii' a=l 

where u^i' change basis from \i) to |q) and uj^ are the density-matrix eigenval- 
ues. The above expression is minimum when using the largest eigenvalues and 
corresponding eigenvectors of the density matrix p, which are all positive or 
zero. So we have shown that the density matrix leads naturally to the optimal 
reduced basis giving the best approximation to the initial state. The last term 
corresponds to the discarded error. 



2.1.1 Quantum information interpretation. Another alternative and very 
interesting interpretation of the density matrix stems from quantum infor- 
mation theory. The DMRG provides an exciting link between strongly cor- 
related systems and quantum information, by providing another perspective 
into quantum phase transitions and interacting wave functions. As stated by 
J. Preskill in j20j "The most challenging and interesting problems in quantum 
dynamics involve understanding the behaviour of strongly-coupled many-body 
systems [...]. Better ways of characterizing the features of many particle entan- 
glement may lead to new and more effective methods for understanding the 
dynamical behaviour of complex quantum systems" 

In Ref. 121] an interpretation of the correlation functions of systems at criti- 
cality is given in terms of wave function entanglement. Recent work 22,23,24^ 
has focused on the quantum information perspective of DMRG (see also |,25j ) . 

To understand the entanglement between the two parts of a bipartite "uni- 
verse", say system S and environment E one can perform a Schmidt decom- 
position of the wave function of the "universe" (assuming it is in a pure 
state), i.e. writing the wave function as a product of states of S and E: 
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min{Ns,NE) 

\i;)= Yl (11) 

a=l 

where Ns and Ne are the Hilbert space sizes of system and environment re- 
spectively. Tracing upon the environment we obtain the system density matrix 

min{Ns,NE) 

PS= w«l«^)(a^l- (12) 

a=l 

Similarly, tracing upon the system degrees of freedom, we obtain pE- Both 
density matrices have the same rank given by r < min{Ns,NE)- One can now 
affirm that the state {tp) is entangled if and only if the Schmidt rank r > 1. 

A quantitative measure of entanglement is given by the von Neumann en- 
tropy |271l28j defined in either subsystem, (here for example in the subsystem 
S), by: 

Ss = -Trps log2 ps = - YZ^^ ^°S2 (13) 

a 

By keeping the highest eigenvalues uj^ of the density matrix, for the general 
case, one obtains the largest entropy S and, hence, the maximum entanglement 
between system and environment. It follows from the singular value decom- 
position theorem that, for a pure target state and any length of system S 
and environment E blocks, Ss + Se + I = Suniverse{= for a pure target or 
universe state), where / is the mutual information of the blocks and measures 
the correlation between them. If both blocks are uncorrelated, then / = 0. It 
also stems from the above that for periodic boundary conditions, where the 
density matrix eigenvalues Wq, have a slower decay, the block entropy is also 
larger than for open boundaries. Therefore, the DMRG will perform better for 
models or setups which have lower block entropy. 

This quantum information perspective leads to an interesting analysis of 
the performance of the DMRG in different dimensions. The entropy S is pro- 
portional to the number of states to be kept in order to maintain the relevant 
information and it depends on system size N (Sn)- Using geometric arguments 
in a (d-l-l)-dimensional field theory including a (d-l)-dimensional hypersurface 
dividing the universe in two (S+E), it is shown that the entropy, which resides 
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essentially at the surface, scales as the hypersurface area PH] 

Sn a {N/Xf-^ (14) 

were A is an ultraviolet cutoff. For one dimension (d = 1) a more detailed 
calculation leads to a logarithmic scaling for gapless critical systems having 
the universal form Sn = c/31n(A^) + Ai (where c is the central charge of the 
underlying conformal field theory and Ai is related to the ultraviolet cutoff). 
One also obtains a saturated entropy for non critical, gapped systems |l24l3UI31j 
when the system size exceeds the correlation length (which is infinite for critical 
systems). This shows that for two dimensions, as the number of relevant states 
to be kept increases with system size, the DMRG isn't an appropriate method 
as it is conceived. It is, however, quite reliable for sufficiently small systems. 

The influence of open boundary conditions on the entropy of entanglement 
for critical XYZ spin chains was analyzed in [22]) finding an additional alter- 
nating term connected with the antiferromagnetic nature of the model. 

There have been several recent publications exploiting the quantum infor- 
mation perspective of the DMRG. For example, the time dependence of the 
the block entropy in spin models with sudden changes in the anisotropies 
was studied in [221 ■ measure the entanglement between two halves of an 
anisotropic Heisenberg chain separated by an impurity bond, the dependence 
of the block entropies (calculated using the density matrix as explained above) 
with system size and anisotropy was done in finding different behaviours 
between the antiferro and ferromagnetic cases. The entropy of two sites in a 
one-dimensional system S = —Tr pj^j+i log2 Pi,i+i was used as a detector of 
phase transitions in spin and fermionic systems and applied to models like the 
bilinear-biquadratic S = 1 Heisenberg chain, the ionic Hubbard model and the 
neutral- ionic transition in a donor-acceptor model for molecular chains |35j . 
A similar quantity was used in 0^ to calculate the entanglement between 
molecular orbitals. This two-site method turns out to work better than the 
single-site entropy for the detection of quantum phase transitions (QPT) pro- 
posed in [H7j since it comprises non local correlations and can detect phases 
with off-diagonal long range order. An even better approach seems to be the 
detection of QPT by using block entropies as proposed in ■38] (without resort- 
ing to the DMRG) as also shown in |37j . 

2.2 General remarks 

With respect to accuracy and convergence, the DMRG behaves in different 
ways depending on the nature of the problem. Its success strongly relies on 
the existence of the so-called matrix product states MPS, a simple version of 
which is the AKLT (AfHeck-Kennedy-Lieb-Tasaki) state for spin 1 chains ,,39] . 
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The DMRG can be viewed as a variational approach within this formula- 
tion |40U41lOn] . For cases where the quantum ground state can be repre- 
sented exactly as a matrix-product state, Ostlund and Rommer showed that 
the renormalization eventually converges to a fixed point and the density ma- 
trix has a finite spectrum of non vanishing eigenvalues (see also |42j . where 
a nice example is given based on the non-hermitean g-symmetric Heisenberg 
model using corner transfer matrices). Gapful models also lead to an excel- 
lent performance of the method, since the density matrix eigenvalue spectrum 
decays exponentially and so does the truncated weight with increasing the 
number of states kept m |44|l4,'-{j. Instead, in the case of critical ID mod- 
els with algebraically decaying correlation functions, it has been shown that 
the eigenvalue spectrum decays much slower |45j which slows down further 
with system size j4()j . An analysis of how to circumvent the problems aris- 
ing near critical points by exploiting conformal field theory predictions using 
multi-targeted sates is done in and, in |3H|) a mathematical description 
of how matrix-product states approximate the exact ground state is given. 
For two-dimensional systems it can be shown that the eigenvalue decay of the 
density matrix is extremely slow, thus the method becomes very unreliable 
for dimensions higher than one [32] • An interesting optimization study in mo- 
mentum space using concepts of quantum information entropy to select the 
best state ordering configurations was performed in Ref. [2Sj and also applied 
to quantum chemistry calculations. They conclude that the optimum results 
are obtained when the states with maximum entropy are located at the center 
of the chain. However, a subsequent study by the same authors led to the 
conclusion that a criterion based on the entropy profile alone does not lead to 
optimal ordering [SHI and further considerations are given based on quantum 
data compression concepts. In Ref. j23j the DMRG's low performance with 
periodic boundary conditions is analyzed and an alternative proposal based 
on highly entangled states is presented which leads to a much higher precision, 
comparable to the open boundary conditions results. In a recent work .51. , S. 
White proposes an interesting correction to the DMRG method which con- 
sists of including only a single site added to the system, where the subsequent 
incompleteness of the environment block is compensated by considering a cor- 
rected density matrix which takes into account fluctuations of the added site. 
Comparing this method to the standard two-added-site method he obtains a 
similar accuracy for open boundary conditions (but faster by a factor of 3) 
and a much higher accuracy for periodic boundary conditions. 

An analytical formulation combining the block renormalization group with 
variational and Fokker-Planck methods is detailed in 52 . The connection of 
the method with quantum groups and conformal field theory is treated in 
|53j and an interesting derivation of the reduced density matrix for integrable 
fermionic and bosonic lattice models from correlation functions was done in 
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|54' . The articles mentioned above give a deep insight into the essence of the 
DMRG method. 

2.3 Symmetries 

It is crucial to include symmetry in the DMRG algorithm since important 
reductions in the Hilbert space size can be achieved by fixing the quantum 
numbers associated with each symmetry. It allows to consider a smaller num- 
ber of states, enhance precision and obtain eigenstates with definite quantum 
numbers. 

For example, most of the models treated with DMRG conserve the total z 
spin projection Sz and particle number in the whole system, so special care has 
to be taken when constructing the superblock Hamiltonian, remembering that 
each block is constituted by several states with different quantum numbers. It 
is easy to see that when the superblock quantum number can be obtained as 
an addition of quantum numbers of each constituent part of the system, the 
density matrix is block diagonal. 

Non abelian or non additive symmetries like SU(2) or total spin are more 
difficult to implement. Sierra and Nishino showed that Interaction-round- a- 
face (IRF) Hamiltonians represent Hamiltonians with a continuous symmetry 
and developed a DMRG algorithm appropriate to this model j55| . This method 
was applied to several spin models with excellent performance |551I56| . Total 
spin conservation, continuous symmetries and parity have been treated also 
in [S7t>60..61..58t„59) and recently used in Ref. in carbon nanotubes, where 
flat band ferromagnetism was found. 

Contrary to momentum-space DMRG (see Sect 5.1), in real space the total 
momentum cannot be fixed, because there is no way to fix the corresponding 
phases in the states that diagonalize the density matrix. However, a very good 
approximation to having states with fixed momenta can be done by keeping 
the appropriate target states [HH] . 

3 Applications 

Since its development in 1992, the number of papers using DMRG has grown 
enormously and other improvements to the method have been performed. For 
example, since 1998 there have been around 80 papers a year using DMRG. 
There have also been several improvements to the method and it is now used 
in areas that are very different to the original strongly correlated electron 
system field. We would like to mention some applications where this method 
has proved to be useful. Other applications related to further developments of 
the DMRG will be mentioned in subsequent sections. 
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An important benchmark in calculations using DMRG was achieved by 
White and Huse 64 when calculating the spin gap in a S* = 1 Heisenberg 
chain obtaining A = 0.41050J with unprecedented accuracy. They also calcu- 
lated very accurate spin correlation functions and excitation energies for one 
and several magnon states and performed a very detailed analysis of the exci- 
tations for different momenta. They obtained a spin correlation length of 6.03 
lattice spacings. Simultaneously S0rensen and Affleck jHS] also calculated the 
structure factor and spin gap for this system up to length 100 with very high 
accuracy, comparing their results with the nonlinear a model. In a subsequent 
paper |66j they applied the DMRG to the anisotropic 5 = 1 chain, obtaining 
precise values for the Haldane gap. They also performed a detailed study of 
the 5 = 1/2 end excitations in an open chain. Thermodynamical properties 
in open 5 = 1 chains such as specific heat, electron paramagnetic resonance 
(EPR) and magnetic susceptibility calculated using DMRG gave an excel- 
lent fit to experimental data, confirming the existence of free spins 1/2 at the 
boundaries ^67^68,69,, 70) . A related problem, i.e. the effect of non-magnetic im- 
purities in spin systems (dimerized, ladders and 2D) was studied in [7T1I72] . In 
addition, the study of magnon interactions and magnetization of 5 = 1 chains 
was done in j73j , supersymmetric spin chains modelling plateau transitions in 
the integer quantum Hall effect in [J^ and ESR studies in these systems was 
considered in [J^. For larger integer spins there have also been some studies. 
Nishiyama and coworkers ^76! calculated the low energy spectrum and corre- 
lation functions of the 5 = 2 antiferromagnetic Heisenberg open chain. They 
found 5 = 1 end excitations (in agreement with the Valence Bond Theory). 
Edge excitations for other values of 5 have been studied in Ref. [77j. Almost 
at the same time Schollwock and Jolicoeur j78j calculated the spin gap in the 
same system, up to 350 sites, (A = 0.085 J), correlation functions that showed 
topological order and obtained a spin correlation length of 49 lattice spacings. 
More recent accurate studies of 5 = 2 chains are found in 79,"5H','SIT and of 
5 = 1 chains in staggered magnetic fields |81i including a detailed comparison 
to the non-linear sigma model in [H^j. In Ref. |H3| the dispersion of the single 
magnon band and other properties of the 5 = 2 antiferromagnetic Heisenberg 
chains were calculated and the phase diagram of S=l bosons in ID lattices, 
relevant to recent experiments in optical lattices was studied in j84j . 

Concerning 5 = 1/2 systems, DMRG has been crucial for obtaining the loga- 
rithmic corrections to the 1/r dependence of the spin-spin correlation functions 
in the isotropic Heisenberg model |HS1- For this, very accurate values for the 
energy and correlation functions were needed. For = 100 sites an error of 
10~^ was achieved keeping m = 150 states per block, comparing with the exact 
finite-size Bethe Ansatz results. For this model it was found that the data for 
the correlation function has a very accurate scaling behaviour and advantage 
of this was taken to obtain the logarithmic corrections in the thermodynamic 
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limit. Other calculations of the spin correlations have been performed for the 
isotropic |86»»87, and anisotropic cases ^88^. Luttinger liquid behaviour with 
magnetic fields have been studied in [HH, field-induced gaps and string order 
parameters in [30] and |3J respectively, anisotropic systems in 15^1^3] and 
the Heisenberg model with a weak link in |^. An analysis of quantum crit- 
ical points and critical behaviour in spin chains by combining DMRG with 
finite-size scaling was done in |25|. Spin systems in more complex lattices like 
the Bethe lattice were considered in jUHl and the effect of twisted boundary 
magnetic fields in 97'j- 

Similar calculations have been performed for the S = 3/2 Heisenberg 
chain jHE]- In this case a stronger logarithmic correction to the spin corre- 
lation function was found. For this model there was interest in obtaining the 
central charge c to elucidate whether this model corresponds to the same uni- 
versality class as the S = \/2 case, where the central charge can be obtained 
from the finite-size scaling of the energy. Although there have been previous 
attempts these calculations presented difficulties since they involved also 
a term ~ 1/ In^ N. With the DMRG the value c = 1 was clearly obtained. 

In Ref. 100 , DMRG was applied to an effective spin Hamiltonian obtained 
from an SU(4) spin-orbit critical state in ID. Other applications were done to 
enlarged symmetry cases with SU(4) symmetry in order to study coherence in 
arrays of quantum dots |101j , to obtain the phase diagram for ID spin orbital 
models |102j and dynamical properties in a magnetic field |103j . 

Dimerization and frustration have been considered in Refs. |1041I10^I1061 
[WlirnHllTIMIllOllllllITT^ and alternating spin chains in [Tn] . 

Several coupled spin chains (ladder models) have been investigated in jll4l 
lll51lll6l lll7''118l lll9U120j . spin ladders with cyclic four-spin exchanges in 
jl21lll22l ll23, 124j and Kagome antiferromagnets in |125j . Zigzag spin chains 
have been considered in jl26lll27lll28j . spin chains of coupled triangles in jl29l 
I1301ll31j . triangular Ising models in jl32j and three-leg spin tubes in |133j . As 
the DMRG's performance is optimal in open systems, an interesting analysis 
of the boundary effect on correlation functions is done in ^ISj . Magnetization 
properties and plateaus for quantum spin chains and ladders |134II135II136II1371 
I138j have also been studied. An interesting review on the applications to some 
exact and analytical techniques for quantum magnetism in low dimension, 
including DMRG, is presented in jl39j . 

There has been a great amount of applications to fermionic systems such 
as ID Hubbard, ionic Hubbard and t-J models [T?ni[TTnrr^l1 4:^111 44111^ 
IIinilIlZlIl48,149,MaiinillE21, Luttinger liquids with boundaries [ESI, the 
Falicov-Kimball model |154j . the quasiperiodic Aubry- Andre chain jl55j and 
Fibonacci-Hubbard models |156|. A recent calculation obtains accurate values 
of the Luttinger-liquid parameter Kp in Hubbard and spin chains from static 
correlation functions |157j . It has also been applied to field theory |lipi58j . The 
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method has been very successful for several band Hubbard models |159[I16U] 
and extended Hubbard models modelling chains of CuO^ plaquettes jl61j . 
Hubbard ladders jMUMlIinilMlIISn] and t-J ladders HHZI- Recent calcu- 
lations in doped Hubbard ladders include the observation of stripes in the 
6-leg ladder for large values of the interaction parameter and not for low val- 
ues |168j , the study of orbital currents and charge density waves in ladders of 
up to 200 rungs jl69j and existence of charge order induced by electron-lattice 
interaction in coupled ladders |170j . Ring exchange on doped two-leg ladders 
were investigated in |171j . Also several coupled chains at different dopings 
have been considered jl721ll73j as well as flux phases in these systems jl74| . 
Time reversal symmetry-broken fermionic ladders have been studied in jl75j 
and power laws in spinless-fermion ladders in jl76j . Long-range Coulomb in- 
teractions in the ID electron gas and the formation of a Wigner crystal was 
studied in jl77j . Several phases including the Wigner crystal, incompressible 
and compressible liquid states, stripe and pairing phases, have been found 
using DMRG for 2D electrons in high magnetic fields considering different 
Landau levels jl781ll79lll80j . Persistent currents in mesoscopic systems have 
been considered in |181j . 

Impurity problems have been studied for example in one- jl82j and two- 
impurity jl83j Kondo systems, in spin chains |184j and in Luttinger Liquids 
and Hubbard models jl851ll86] . There have also been applications to Kondo 
and Anderson lattices [TK7lllH8lll?^<)lll^()lll< jl . l9lj l<j:-il[Tnil[^llf)Klllf)6j . Kondo 
lattices with localized configurations [l97j . the two-channel Kondo lattice 
on a ladder jl98j and on a chain (199] , a t-J chain coupled to localized Kondo 
spins '20Qj and ferromagnetic Kondo models for manganites |2()ip2()2p2()3p2()4j . 
See also Refs. I205U206] for recent results for the antiferro and ferromagnetic 
Kondo lattice models. 

The DMRG has also been generalized to ID random and disordered systems, 
and applied to the random antiferromagnetic and ferromagnetic Heisenberg 
chains |207|l208j . including quasiperiodic exchange modulation (20^j and a 
detailed study of the Haldane 210 and Griffiths phase |211j in these systems. 
Strongly disordered spin ladders have been considered in |212j . It has also 
been used in disordered Fermi systems such as the spinless model |213[I214] . 
In particular, the transition from the Fermi glass to the Mott insulator and 
the strong enhancement of persistent currents in the transition was studied in 
correlated one-dimensional disordered rings |215j . Disorder-induced crossover 
effects at quantum critical points were studied in j21fij- 

There have been recent applications of the DMRG in nanoscopic devices 
such as transport properties through quantum dots j217j . the study of the 
influence of interactions in a reservoir on the levels of a quantum dot cou- 
pled to it j218j and charge sensing in quantum dots ^219) . A very interesting 
application of the DMRG to study single wall carbon nanotubes, where the 
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tube is mapped onto a ID chain with longer range interactions (depending 
on the quirahty) can be found in Ref. (220^ and the study of edge states 
in doped nanostructures in |221j . Inspired in recent experiments of trapped 
bosonic atoms |222j . a one dimensional Bose gas was studied using ab-initio 
stochastic simulations at finite temperatures covering the whole range from 
weak to strong interactions |223j . They used a discretized version of the model 
and a block factorization of the kinetic energy, using the density matrix to 
select the most relevant states from each block. Also models for trapped ions 
have been studied using DMRG, where different quantum phases were ob- 
tained |224j ■ A recent study concerns the colossal electroresistance (due to 
band bending) through the interface between a metal and a strongly interact- 
ing system (in ID) in the presence of an electric potential obeying the Poisson 
equation. In this case the use of the finite-size method played an important 
role in the convergence of the ground state due to the local density-dependent 
potential. |225j 



4 Bosonic degrees of freedom 

A vast research area in correlated systems concerns bosonic degrees of freedom. 
For example, phonons play fundamental role in models of currently interest- 
ing correlated systems such as high-Tc superconductors, manganites, organic 
compounds, or nanoscopic systems. 

A significant limitation to the DMRG method is that it requires a finite basis 
and calculations in problems with infinite degrees of freedom per site require 
a large truncation of the basis states |226j . However, Jeckelmann and White 
developed a way of including phonons in DMRG calculations by transforming 
each boson site into several artificial interacting two-state pseudo-sites and 
then applying DMRG to this interacting system |227j (called the "pseudo-site 
system"). The proposal is based on the fact that DMRG is much better able to 
handle several few-states sites than few many-state sites |228j . The key idea is 
to substitute each boson site with 2^ states into N pseudo-sites with 2 states 
|229j . They applied this method to the Holstein model for several hundred sites 
(keeping more than a hundred states per phonon mode) obtaining negligible 
error. In addition, to date, this method is the most accurate one to determine 
the ground state energy of the polaron problem (Holstein model with a single 
electron). For a recent comprehensive review on exact numerical methods in 
electron-phonon problems see j23Uj . 

This method has been applied recently to the calculation of pairing, 
CDW and SDW correlations in the Holstein-Hubbard model with competing 
electron-electron and electron-phonon interactions, finding an enhanced su- 
perconducting exponent when next-nearest neighbour hopping is added |231j . 
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Figure 5. Ground-state energy of a 4 site, half-filled Holstein system vs. the number of phonon 
states kept on each site of the lattice, for the optimal and bare phonon basis. Reprinted from 12321 

with permission. 

They propose here a way to improve the bare repulsion felt by the electrons in 
the newly added pseudo-sites in the infinite-size algorithm by adding a chemi- 
cal potential to these sites, finding lower ground state energies. An interesting 
result was obtained in |233j using this method combined with the dynamical 
DMRG (see Sect 18. 1(1 to calculate single-particle spectra in one-dimensional 
Mott insulators (the Holstein-Hubbard model at half filling in this case) . They 
observe that the charge-spin separation characteristic of ID systems is robust 
against electron-phonon coupling. 

An alternative method (the "Optimal phonon basis") |232j is a procedure 
for generating a controlled truncation of a large Hilbert space, which allows 
the use of a very small optimal basis without significant loss of accuracy. 
The system here consists of only one site and the environment has several 
sites, both having electronic and phononic degrees of freedom. The density 
matrix is used to trace out the degrees of freedom of the environment and 
extract the most relevant states of the site in question. In following steps, 
more bare phonons are included to the optimal basis obtained in this way. 
This method was successfully applied to study the Holstein model (FigEJ, the 
interactions induced by quantum fluctuations in quantum strings as an ap- 
plication to cuprate stripes |234j . the dissipative two-state system |235j and, 
recently, the spin-Peierls model |236j . A variant of this scheme is the "four 
block method", as described in j237j . They obtain very accurately the Lut- 
tinger liquid-CDW insulator transition in the ID Holstein model for spinless 
fermions. 

The method has also been applied to pure bosonic systems such as the 
disordered bosonic Hubbard model j238j . where gaps, correlation functions 
and superfiuid density are obtained. The phase diagram for the non-disordered 
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Bose-Hubbard model, showing a reentrance of the superfluid phase into the 
insulating phase was calculated in Ref. |239j . It has also been used to study 
a chain of oscillators with optical phonon spectrum j24Uj and optical phonons 
in the quarter-filled Hubbard model for organic conductors [241j . 



5 Generalized DMRG 

Since its development, the DMRG has been widely improved in several direc- 
tions in order to deal with different kinds of problems. It has been generalized 
to treat models where a representation in momentum or energy space is more 
appropriate. Based on these ideas, an important extension towards Quantum 
Chemistry and the calculation of physical properties of small grains and nuclei 
was performed. Below we describe each of these ideas. 

5.1 Momentum representation 

In 1996, Xiang j242j proposed an alternative formulation of the DMRG in 
momentum space, rather than in real space, which can be used also in two or 
more dimensions, it deals better with periodic boundary conditions and has 
the advantage that the total momentum quantum number can be fixed. This 
algorithm has inspired the extension of the DMRG to finite fermionic systems, 
small grains and nuclear physics. The main idea is to write the Hamiltonian 
in momentum space. If the original model is long ranged, as with Coulomb 
interactions, then it will turn local in the reciprocal space and the DMRG 
could be used in the normal way, where now each "site" is a state with definite 
momentum. The most convenient growth procedure for real space short-ranged 
Hamiltonians is to start using one-particle states lying near the Fermi energy, 
above and below, and subsequently add additional more distant states. For 
longer ranger models other orderings are to be considered |28j . 

In his original work, Xiang applied these ideas to the ID and 2D Hubbard 
models with nearest-neighbour interaction, defined in general as: 



The operator cj^ (cjo-) creates (annihilates) an electron on site i with spin 
projection a, riia = cl^Cia, t is the hopping parameter and U the on-site 
Coulomb interaction. 
In momentum space it reads: 




(15) 
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where nko- = cj^o-CkCT is the number operator for particles with spin a and 
momentum k and ek is the energy dispersion of the particles. 

By defining appropriate operator products using the original ones, the num- 
ber of operators to be kept can be greatly reduced and the convergence be- 
comes faster and more accurate. Other authors j243ll28j have improved the 
method and performed more systematic and interesting analysis on the order 
of states to be kept, depending on the kind of interactions in the model and 
very accurate dispersion relations have been obtained in |173j . The method 
has a better performance than the real space DMRG in two dimensions for 
low values of the interaction (as seen in 4x4 lattice calculations), the accuracy 
of the method becoming worse for larger values of U. Recently, an interesting 
algorithm for the calculation of eigenstates with definite momentum has been 
put forward in |244j . based on the PEPS (projected entangled pair states). 
In this case additional degrees of freedom are added for each site and the 
particular momentum quantum number is fixed by constructing a "projected 
entangled-multiparticle state" (PEMS). It has been successfully applied to cal- 
culate excitation spectra in the bilinear-biquadratic 5 = 1 Heisenberg model. 

5.2 Grains, nuclei and high energy physics 

Another successful extension to the DMRG is the calculation of physical prop- 
erties of small fermionic systems, like superconducting grains and models for 
nuclei. This was developed by Dukelsky, Sierra and collaborators |245j . For a 
review on these applications see the recent publication _.246_. 

A typical and general Hamiltonian for these finite systems consists of one- 
and two-particle terms of the form 



where Tij represents the kinetic energy and Vijki the residual two-body inter- 
action between the effective particles. 

Starting from the single particle levels and a certain filling which defines 
the Fermi energy, one can define particle and hole blocks in a similar manner 
as the real space DMRG (see Fig. The empty and occupied levels are the 
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particle levels ■ 



hole levels 



Figure 6. DMRG method for finite fermionic systems where the added "sites" he successively at 

larger distances from the Fermi energy. 



particle and hole levels respectively (referring to the character of the single- 
particle excitations). By successively adding levels at increasing distance from 
the Fermi level, the DMRG renormalization is applied. This might not be the 
appropriate level ordering for DMRG for more complex Hamiltonians |28ll246j . 

Within this scheme, superconductivity in small grains has been studied using 
the half-filled reduced BCS Hamiltonian 
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where i,j = 1,2, label the single-particle energy levels, Cja are the 

electron annihilation operators associated with the two time-reversed states 
|j, C7 = -|-) and |j, o" = — ) , ^ is the chemical potential, A is an adimensional 
BCS coupling constant and d is the single-particle energy level spacing (as- 
sumed to be constant). By using DMRG a smooth superconducting-normal 
state crossover with decreasing system size was found, a more reliable re- 
sult than the abrupt transition obtained using exact diagonalization and BCS 
approximation. This method was also used for two superconducting grains 
coupled by tunneling j247j . 

Other recent applications have been in nuclear shell model calculations 
where a two-level pairing and pairing plus quadrupole interactions model have 
been considered [248] . Here the nucleus is modelled by completely filled core 
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shells and valence orbitals partially filled with protons and neutrons where 
the starting states are the Hartree Fock functions. A realistic calculation for 
the ^^Mg nucleus showed very slow convergence mainly due to the lack of 
improvement of the finite-system algorithm in this basis and the fact that 
total angular momentum cannot be explicitly fixed. A recent work obtained 
the sd- and /p-shell models for ^^Si and ^^Ni respectively in a different basis 
using proton-neutron wave function factorization which allows fixing total an- 
gular momentum plus the finite-system algorithm '249^ The DMRG has also 
been extended to study complex-symmetric density matrices which describe 
many-body open quantum systems (where resonant and non-resonant levels 
are considered) and applied to the unbound nucleus ^He with possible further 
applications to open quantum dots, open microwave resonators, etc. j25Qj 

A very interesting and successful application is a recent work in High Energy 
Physics |251j . Here the DMRG is used in an asymptotically free model with 
bound states, a toy model for quantum chromodynamics, namely the two 
dimensional delta-function potential. For this case an algorithm similar to 
the momentum space DMRG |242j was used where the block and environment 
consist of low and high energy states respectively. The results obtained here are 
much more accurate than with the similarity renormalization group |252j and 
a generalization to field-theoretical models is proposed based on the discreet 
light-cone quantization in momentum space |253j . 

5.3 Molecules and Quantum Chemistry 

Using the standard real space DMRG, there have been several applications to 
molecules and polymers, such as the Pariser-Parr-Pople (PPP) Hamiltonian 
for a cyclic polyene j254|l255] (where long-range interactions are included), 
magnetic Keplerate molecules |25f)| . molecular Iron rings 257^, spin crossover 
molecules |258j and polyacenes considering long range interactions |259j . The 
application to conjugated organic polymers was performed by adapting the 
DMRG to take into account the most important symmetries in order to obtain 
the desired excited states (MI- Also conjugated one-dimensional semiconduc- 
tors |26nj have been studied, in which the standard approach can be extended 
to complex ID oligomers where the fundamental repeat is not just one or two 
atoms, but a complex molecular building block. Relatively new fields of ap- 
plication are the calculation of dynamical properties in the Rubinstein-Duke 
model for reptons j261j and excitons in dendrimer molecules j262j . 

Another very important problem is the ab initio calculation of electronic 
states in Quantum Chemistry. The standard approaches consist of the Den- 
sity Functional Theory and the Hartree Fock starting point to calculate com- 
plex electronic structure (using for example the Configuration Interaction ap- 
proach). When electron interactions are included, relatively small clusters can 
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Figure 7. Ground state energy (left) and error (right) of a water molecule in a 25 orbital basis 
using several methods. The DMRG results are plotted as a funtion of the number of states kept, m 
(not relevant for the other energies). Reprinted from |263| with permission. 



be considered due to the large Hilbert spaces involved. To be able to treat 
larger molecules, in 1999, White and Martin j263j proposed a way of using 
the DMRG for quantum chemistry calculations, based on the momentum rep- 
resentation DMRG |242j . They considered a general model for molecules, Eq 
(|17l) . where the electron- nucleus interaction in considered in the first term 
and the operators indices represent orbitals in the Hartree Fock basis (see 
also |,264.265t,.266) ). Here, DMRG is applied within the conventional quantum 
chemical framework of a finite basis set with non-orthogonal basis functions 
centered on each atom. After the standard Hartree-Fock (HF) calculation in 
which a Hamiltonian is produced within the orthogonal HF basis, DMRG is 
used to include correlations beyond HF, where each orbital is treated as a 
"site" in a ID lattice. One important difference with standard DMRG is that, 
as the interactions are long-ranged, several operators must be kept, making 
the calculation somewhat cumbersome. However, very accurate results have 
been obtained in a check performed in a water molecule (keeping up to 25 
orbitals and m ~ 200 states per block), obtaining an offset of 2.410~^Hartrees 
with respect to the exact ground state energy J267 , a better performance than 
any other approximate method |263j (see FigEJ. 

In order to avoid the non-locality introduced in the treatment explained 
above. White introduced the concept of orthlets, local, orthogonal and com- 
pact wave functions that allow prior knowledge about singularities to be in- 
corporated into the basis and an adequate resolution for the cores |264j . The 
most relevant functions in this basis are chosen via the density matrix. An 
application based on the combination with the momentum version of DMRG 
is used in |265j to calculate the ground state of several molecules. In Ref. (268) . 
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the efficiency of the DMRG method was analyzed. Although it works very well 
for some molecules, like Be2, HF, H2O, CH4, calculations in the N2 molecule 
were not satisfactory. The reasons for the failure could be related to the level 
ordering taken for the DMRG algorithm, although it is not conclusive. A very 
careful analysis of the method presented in a pedagogical way can be found in 
Ref. P69.. 

In a recent work, quantum information theory concepts were used together 
with DMRG to calculate orbital interaction in different molecules, such as LiF, 
CO, and F2 - Here the interaction between electrons in two orbitals in a 
many-body wave function was calculated as the difference between the entan- 
glement of both orbitals taken together and individually (always belonging to 
the same DMRG "system") for a given wave function. An interesting analysis 
of the relevance of the level ordering within the DMRG growing procedure is 
given. 



6 Higher dimensions {D > 1) 

The possibility of using DMRG in two dimensions is currently one of the 
most challenging problems. As we mentioned above equation 114(1 . due to the 
increasing number of states to be kept in dimensions higher than one, the 
DMRG is not suited for these cases. By analyzing analytically the spectra 
of a two-dimensional system of coupled oscillators, Chung and Peschel j49j 
showed that the density matrix eigenvalues decay much slower than in the 
one-dimensional case. However, if the system size is kept relatively small, in 
two dimensions one can reach larger systems than with exact diagonalization 
methods and the zero temperature results are more reliable than with Monte 
Carlo methods. The first attempts to extend the DMRG to two dimensions 
considered different one-dimensional kind of paths like the one shown in figure 
IHl where the DMRG is used in the standard way. The main drawback of this 
method is that many operators have to be kept and the interactions are long- 
ranged in this ID mapping. 

Quite large quasi-2D systems can be reached, for example in |27nj where a 
4x20 lattice was considered to study ferromagnetism in the infinite-U Hubbard 
model; the ground state of a 4- leg t-J ladder in |271j : the one- and two- hole 
ground states in 9x9 and 10x7 t-J lattices in |272j : a doped 3-leg t-J ladder 
in |273j : the study of striped phases in j274j : domain walls in 19x8 t-J systems 
in j275j ; the 2D t-J model in P76l| and the magnetic polaron in a 9x9 t-J lattice 
in |277j . Also big CaV^Og spin-1/2 lattices reaching 24x11 sites have been 
studied as well as a recent study of the spin liquid phase in the anisotropic 
triangular Heisenberg model in up to 8x18 lattices j278j . 

There have also been some recent alternative attempts to implement DMRG 
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Figure 8. General sweeping and superblock configuration for a two-dimensional lattice (from 1131 . 

with permission) 



in two and higher dimensions with clever block arrangements j2791l243ll28Ul 
I2811l282ll283j but the performances are still poorer than in ID. A recent ex- 
tension using a novel DMRG algorithm for highly anisotropic spin systems 
has shown promising results |284j . Here a two-step method is proposed which 
is based on the diagonalization of a single chain to obtain the lowest-lying 
states and then using DMRG to couple these chains in a ID manner. As long 
as the energy width of the retained states is much larger than the interchain 
coupling, this method yiels results which are comparable to Quantum Monte 
Carlo calculations. 

Another problem in physics that has been handled using DMRG is that of 
electrons in a high magnetic field. Shibata and Yoshioka jl78-179f,180) used the 
eigenstates of free electrons in a perpendicular magnetic field in the Landau 
gauge to represent the local orbitals. By choosing these single particle states as 
effective sites, where the wave functions are parametrized by only two quantum 
numbers (the Landau level and the x component of the center of coordinates 
of an electron in a cyclotron motion), it provides a natural mapping onto a 
ID system (see also j285j ). The long-range Coulomb interactions between the 
electrons stabilizes various different electronic states depending on the filling 
of the Landau levels. If the magnetic field is larger than the characteristic 
Coulomb interactions, the electrons in filled Landau levels can be considered 
inert and the ground state is determined by the partially filled level. The 
ground state properties for different filling factors were characterized using 
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the expectation value of the pair correlation function 

^ " iVe ^n/- 1) E + - R^-) 1^) ' (19) 

where Rj is the center coordinate of the ith electron's wave function, r is 
the relative distance between the pair of electrons and L^^y are the unit cell 
dimensions. In jl78j the phase diagram for the third Landau level was obtained, 
showing a Wigner crystal for low filling factors, up to u = 0.24, where a 2- 
electron bubble phase forms. For > 0.38, a striped phase shows up (see Fig. 
EI) and in |18()j an analysis of the excitonic behaviour of the = 1 bilayer 
quantum Hall systems was done. 



(a) Wigner crystal 



1.65- 



(b) two-electron 
bubble 
-12.70 I » ^ 

Hi 



-6.35 - 



.V 



6.35- 



12.70- 



— I — r 

°° ^ o 3: S 

cn y ^ ^ 



(c) Stripe 





Figure 9. Pair correlation function of a) a striped state, b) a two-electron bubble state and c) a 
Wigner crystal (from 285 , with permission) 



A very promising approach towards two dimensional calculations was re- 
cently proposed by Verstraete and Cirac |22j, based on the projected entangled- 
pair-states (PEPS). It relies on the matrix-product states (MPS) approach to 
describe the variational wave functions which is the basis of success of the 
DMRG. Following Ref. |22]'s notation, let's represent the ID state \"^) as a 
MPS of N d-dimensional physical systems \s) (e.g. d = 2 for individual 5 = 1/2 
spins). We then enlarge the Hilbert space by replacing each physical site \s) 
by two auxiliary systems and each of a certain dimension D (to be ob- 
tained variationally to minimize the energy), except for the extremes of the 
chain, where only one auxiliary system is added, and set them in the highest 
entangled state |(/)) = 'Ylin=i (see Fig. [TU)) . Each of these states is then 
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Figure 10. Schematic representation of MPS (a) 1 dimension and (b) 2 dimensions, (c) 
representation of 2D PEPS. The bonds represent pairs of maximally entangled D— dimensional 
auxiliary spins and the circles denote projectors that map the auxiliary spins to the physical ones. 

Reproduced from fZZ] with permission. 



projected to the real spin by applying an operator Qfc to each pair to obtain 
the real entangled system, 



1^) = Qi(^Q2® ...Qn \(t>)---\4>) 

d 

= Yl F,{Al\...,A'^n\si,...,SN). (20) 

Sl,...,SiV = l 

Here the matrices Af. have elements (A^)/^^ = {slQkll,''') and the indices / and 
r of each matrix Af. are related to the left and right bonds of the auxiliary 
systems with their neighbors. The function Fi denotes the trace of the product 
of the matrices. An extension of this state was proposed for two dimensions 
(see Fig. llDh). where now each site is represented as the projection of 4 D- 
dimensional auxiliary systems. The ground state of a Hamiltonian corresponds 
to the PEPS with a given dimension D which minimizes the energy and this 
is performed by iteratively optimizing each A tensor while keeping the others 
fixed, in a similar manner as the finite-size DMRG zipping procedure. Within 
this MPS picture, it was shown also by Verstraete and collaborators j286j that 
a uniform picture can be given of Wilson's numerical renormalization group 
method (8, together with White's DMRG. 
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7 Time-dependent analysis 

One of the most important recent developments of the DMRG concerns real- 
time calculations: the so-called time-dependent DMRG. Its potential applica- 
tions include the analysis of the evolution of wave functions, the calculation 
of transport properties in low-dimensional systems, non-equilibrium transport 
through nanosystems or dissipative quantum mechanics. A recent comprehen- 
sive review on the time-dependent DMRG can be found in Ref. j287j . 

There are two main different schemes to tackle time-dependent problems. 
Both different approaches have been classified as dynamic (Hilbert space) and 
adaptive time-dependent methods. In the first case the Hilbert space for the 
initial wave function is kept as large as possible so as to be able to describe the 
time-evolved state for sufficiently long times. In the second approach, instead, 
the states kept are calculated at each time interval, adapting to the evolving 
Hamiltonian. 

The first method was developed by Cazalilla and Marston j288j in a sim- 
plified version which targeted only the ground state (defined as static method 
in j287j ). After growing the system up to a certain size, a quantum dot or 
a junction is added and a time-dependent perturbation is H'{t) set on. The 
time-dependent Schrodinger equation is then integrated forward in time using 
a discretized algorithm: 



i^^l^'(i)) = [(//trunc - ^O) + H'{t)]\^{t)) . (21) 

where the initial state is chosen to be the ground state of the unperturbed 
system. As this method does not conserve unitarity, the alternative Crank- 
Nicholson procedure (289) . 

mt + dt)) ^ l-^-^t.unc(tMt/2 
1 + tHtrunc[t)0t/2 

can be used. The calculation of the denominator can be calculated using the 
biconjugate gradient approach |29nj . If more than one target state is used 
during the growing procedure, a higher precision is obtained and the time 
evolution can be followed longer |291j . Fig. 1111 shows the current through two 
quantum systems calculated using this method. 

A slightly different approach was used by Schmitteckert j292j . Instead of us- 
ing the differential equation, he applied the time evolution operator exp(—iHt) 
directly by using a Krylov basis expansion of the matrix exponential, using 
a Lanczos procedure, within the finite-system algorithm. Using this method, 
and keeping information of the evolved wave function in a similar manner as 
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Figure 11. Current through a) a non-interacting quantum dot and b) a junction between two leads 
(defined in 12881 '). Nt is the number of excited target states, L the number of sites and M the 
number of states kept (reprinted with permission from 291 ). 

described below, transport properties of spinless electrons through strongly 
interacting ID systems with arbitrary bias was calculated |293j . This method 
does not involve Trotter approximations (see below) and is appropriate to 
deal with general geometries and interactions. An improvement to this dy- 
namic method was performed in |294j by using the 4th order Runge-Kutta 
method and targeting a small interval of time (time-step targeted method) . 

The adaptive approach |295ll29Uj was based on the algorithm developed 
by Vidal for the time evolution of weakly entangled quantum states |296j . 
originally formulated in the matrix product language. 

For a one dimensional system with nearest-neighbour interactions, the in- 
finitesimal time evolution operator exp{—iHdt) can be decomposed using a 
second order (of course higher orders can also bee applied) Trotter-Suzuki 
expansion into link operators (to order dt^) as: 

^-idtH ^ ^-idtHt/2^-idtH2/2 ^-idtH2 /2 ^-idtHi /2 ^23) 

where Hi is the Hamiltonian term acting on link i. When applied t/dt times 
to get one time interval, the error is of the order {dt^)t, so, in the worst case 
it scales linearly with t (an extensive error analysis is performed in j297j ^. 
One important characteristic of the DMRG is that, when operators are cal- 
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culated acting on the single sites outside the renormalized blocks, the result is 
exact. Taking advantage of this, the above decomposition can be applied se- 
quentially on these exact sites, in a finite-sweep-like procedure. This algorithm 
can be easily incorporated into a static DMRG code. Using this method. White 
and Feiguin j295j calculated the same transport properties as in Fig. ^2 and 
their results match the exact result (full line) very precisely for larger times. 
They have also calculated other time-dependent responses like the dispersion 
of a spin excitation in a 200-sites spin-1 Heisenberg chain (see Fig. I12|) . 





, — — 









Figure 12. Time evolution of a spin excitation wave-packet in a spin-1 Heisenberg chain with 200 
sites. In b) only the 50 leftmost sites are shown. Reprinted from |295| with permission. 



The results obtained by this method are very accurate thus allowing to 
obtain frequency-dependent quantities by Fourier transformation j295,. Some 
recent applications include the calculation of the evolution of one-dimensional 
density waves |298j and the dynamics of the superfluid-to-Mott transition |299j 
(here using the matrix formalism j296j ) in ultracold bosons in an optical lattice 
and the analysis of charge spin separation in cold Fermi gases (300] by observ- 
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ing the time evolution of wave packets |3Ulj . It was also applied to calculate 
the Zener breakdown in band and Mott insulators |3U2j . real-time dynamics 
of spin 1/2 Heisenberg chains |297j . conductance through strongly correlated 
systems |303j and a single atom transistor in ID optical lattice |304j . An in- 
teresting application of this time-dependent analysis is the recent study of the 
evolution of the entropy of a block of spins in the Heisenberg chain after a sud- 
den quench in the anisotropy |S1 , which gives information of the propagation 
of entanglement and critical behaviour. 

This method is very precise and easy to implement. However, it has the dis- 
advantage of being restricted by geometry and nearest-neighbour interactions. 
For other geometries the dynamic method mentioned above is preferred. 



8 Zero temperature dynamics 

DMRG can also be used to calculate dynamical properties of low-dimensional 
systems, useful to interpret experimental results from, for example, nuclear 
magnetic resonance (NMR), neutron scattering, optical absorption and pho- 
toemission, among others. There have been two main approaches to the dy- 
namics: the Lanczos [BHll3U5j and correction vector techniques j3U51l3U6ll3U7j . 
The first gives complete information of the whole excitation spectrum at the 
expense of less accuracy for large systems, specially at high energies. The 
latter, instead, focuses on particular energy values and gives more precise in- 
formation, being numerically much more expensive. 



8.1 Lanczos technique 

We want to calculate the following dynamical correlation function at T = 0: 

CA{t-t') = {MA\t)A{t')\,Po), (24) 

where A'^ is the Hermitean conjugate of the operator A, A{t) is the Heisenberg 
representation of A, and |V'o) is the ground state of the system. Its Fourier 
transform is: 

Ca{oj) = KVnIAlV'o)!' <5(a; - {En - Eo)), (25) 

n 

where the summation is taken over all the eigenstates \^n) of the Hamiltonian 
H with energy £'„, and £^0 is the ground state energy. 
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Defining the Green's function 

GAiz) = {MAHz-H)-'A\ijo), (26) 
tlie correlation function C^(u;) can be obtained as 

Ca{u;) = -- lim Im Ga{uj + irj + Eq). (27) 

Tlie function Ga can be written in the form of a continued fraction: 

GA^>= '^'i-^'-^'t;' (28) 

z-ao 

z—ai — ^_ 

The coefficients a„ and hn can be obtained using the following recursion equa- 
tions ,308.^: 

|/„+i) = - an\!n) - bl\fn~i) (29) 

where 

l/o)=^|V'o> 

an = {fn\H\fn)/{fn\fn), 

bl = {fn\fn)/{fn-l\fn-l); ^0 = (30) 

For finite systems the Green's function (j'a(-z) has a finite number of poles. 
Within this formulation it is not necassary to compute all coefficients a„ and 
bn since the spectrum is nearly converged already with the first few tens or 
hundreds of them, depending on the problem. The DMRG technique presents 
a good framework to calculate such quantities. With it, the ground state, 
Hamiltonian and the operator A required for the evaluation of Ca{(^) are 
obtained. An important requirement is that the reduced Hilbert space should 
also describe with great precision the relevant excited states \ipn)- This is 
achieved by choosing the appropriate target states. For most systems it is 
enough to consider as target states the ground state iV'o) and the first few \fn) 
with n = 0, 1... and |/o) = ^iV'o) as described above. In doing so, states in the 
reduced Hilbert space relevant to the excited states connected to the ground 
state via the operator of interest A are included. The fact that |/o) is an 
excellent trial state, in particular, for the lowest triplet excitations of the two- 
dimensional antiferromagnet was shown in Ref. |31L)j . Of course, if the number 
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Figure 13. Lower bound of the continuous excitation band of the AF 5 = 1/2 Heisenberg model 
for different system sizes L (reprinted from Ref. 305 with permission). The full line is the 
analytical result for the infinite system using Bethe Ansatz. 



m of states kept per block is fixed, the more target states considered, the less 
precisely each one of them is described. An optimal number of target states 
and m have to be found for each case. Due to this reduction, the algorithm 
can be applied up to certain lengths, depending on the states involved. For 
longer chains, the higher energy excitations will become inaccurate. Proper 
sum rules have to be calculated to determine the errors in each case. 

This method has been successfully applied to a number of problems, like spin 
1/2 (see FigCni, 3/2 [HElinHlEHlffl and spin 1 chains OT, the spin-boson 
model '311j, the Holstein model |3121I233] and spin-orbital chains in external 
fields [313) . It was also applied to extract spin-chain dispersion relations |314j . 
dynamics of spin ladders |123j . spectral functions in the infinite-U Hubbard 
model |141j , optical response and spinon and holon excitations in ID Mott in- 
sulators j315ll316] . In a recent work, this method, together with the correction 
vector described below, was used to calculate the spectral function of ultra- 
small Kondo systems and analyzed finite-size and even-odd effects [31 7| with 
great accuracy. In Section 10 we will describe its application as the impurity 
solver within the Dynamical Mean Field Theory method (DMFT). 



New Trends in Density Matrix Renormalization 



35 



8.2 Correction-vector method 

This method focuses on a particular energy or energy window, allowing for a 
more precise description in that range and the possibility of calculating spectra 
for higher energies. Instead of using the tridiagonalization of the Hamiltonian, 
but in a similar spirit regarding the important target states to be kept, the 
spectrum can be calculated for a given z = w + ir] hy using a correction vector 
(related to the operator A that can depend on momentum q). 

Following the Green's function given above, the (complex) correction vector 
\x{z)) can be defined as: 

\x{z)) = ^^^m (31) 

so the Green's function can be calculated as G{z) = {il)o\A'f\x{z)). 

Separating the correction vector in real and imaginary parts \x{z)) = 
\x'^{z)) + i\x'^{z)) we obtain 

((H - wf + 7]^)\x\z)) = -r^Alijo) 

Ix'iz)) = -{w-H)\x'iz)) (32) 
V 

The former equation is solved using the conjugate gradient method. In order to 
keep the information of the excitations at this particular energy the following 
states are targeted in the DMRG iterations: The ground state \ipo), the first 
Lanczos vector ^IV'o) and the correction vector \x{z)). Even though only a 
certain energy is focused on, DMRG gives the correct excitations for an energy 
range surrounding this particular point so that by running several times for 
nearby frequencies, an approximate spectrum can be obtained for a wider 
region 13051. 

A variational formulation of the correction vector technique has been devel- 
oped in |3fl7j . In order to solve Eq. the following equation is minimized 
w.r.t \X): 

WA,r,{uJ, X) = {X\{H - wf + rf\X) + T^iMMX) + r]{X\A\^o). (33) 

For any r] ^ and finite uj this function has a well defined minimum for the 
quantum state which is solution of H32|). i.e. \x^{z)). This method is very similar 
to thge correction vector one, but has a smaller error in the determination of 
the spectral function |23Uj . It has been successfully applied to calculate the 
optical conductivity of Mott insulators '318'. l319ll32Uj . spectral functions of 
the ID Hubbard insulator |316j (see Fig. 114(1 and away from half-filling (321j 
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Figure 14. Single-particle spectral function for the ID Hubbard insulator for U = H and several 
momenta. The black curves are the DMRG results for T = (64 sites) and the blue curve 
corresponds to Quantum Monte Carlo data for T = 0.25f (reprinted from Ref. i316i with 

permission). 



and, more recently, quantum impurity spectral functions for one |,S22j and two 
impurities |323j . 

The correction vector model has also been applied to determine the nonlinear 
optical coefficients of Hubbard chains and related models |324j . to calculate ac 
conductivity of the Bose-Hubbard model 325 and the single-impurity Ander- 
son model. |326j It has also been recently applied to calculate the conductance 
through an interacting spinless fermionic systems in the linear response using 
the Kubo formula |327j . 

An interesting development for calculating response functions in single im- 
purity systems in the presence of a magnetic field was done in |328| by using 
the DMRG within Wilson's NRG to obtain the Green's function. Very recently, 
Verstraete and collaborators proposed a variational method (by optimizing the 
correction vector) for evaluating Green's functions applicable to ID quantum 
chain models, based on the matrix product states formalism and successfully 
applied to the single impurity Anderson model ^286) . 
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9 Classical systems 

The DMRG has been very successfuhy extended to study classical systems. 
For a detailed description we refer the reader to Ref. |329j . Since ID quantum 
systems are related to 2D classical systems j33U| , it is natural to adapt DMRG 
to the classical 2D case. This method is based on the renormalization group 
transformation for the transfer matrix T (TMRG) |331j . It is a variational 
method that maximizes the partition function using a limited number of de- 
grees of freedom, where the variational state is written as a product of local 
matrices [lUI. For 2D classical systems, this algorithm is superior to the clas- 
sical Monte Carlo method in accuracy, speed and in the possibility of treating 
much larger systems. A recent improvement of this method considering peri- 
odic boundary conditions is given in [532_ and a detailed comparison between 
symmetric and asymmetric targeting is done in f,333_. TMRG has also been 
successfully used to renormalize stochastic transfer matrices in a study of cel- 
lular automatons |334j . The calculation of thermodynamical properties of 3D 
classical statistical systems has been proposed |28Uj where the eigenstate of 
the transfer matrix with maximum eigenvalue is represented by the product 
of local tensors optimized using DMRG. 

A further improvement to this method is based on the corner transfer ma- 
trix IS3S1, the CTMRG [SMlESIlESHllSSSEinilMIl and can be generahzed to 
any dimension 1342' . In Ref. |343j . there is an interesting analysis on the low 
effective theory underlying the DMRG, based on a real-space RG of the corner 
Hamiltonian and applied to spin chains. 

It was first apphed to the Ising ^iM^M^i^M and Potts models [SZl, 
where very accurate density profiles and critical indices were calculated. Fur- 
ther applications have included non-hermitian problems in equilibrium and 
non-equilibrium physics. In the first case, transfer matrices may be non- 
hermitian and several situations have been considered: a model for the Quan- 
tum Hall effect |348j . the g-symmetric Heisenberg chain related to the con- 
formal series of critical models |349) and the anisotropic triangular nearest 
and next-nearest neighbour Ising models |131j . In the second case, the adap- 
tation of the DMRG to non-equilibrium physics like the asymmetric exclusion 
problem ^_350_ and reaction-diffusion problems j42ii351i| has shown to be very 
successful. It has also been applied to stochastic lattice models like in |352j and 
to the 2D XY model |353j . to the study of simplified models for polimerized 
membranes in thermal equilibrium (folding of triangular lattices) |354| . 
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10 Finite-temperature DMRG 

The first attempt to include the effect of finite temperature in the DMRG 
procedure was performed by Moukouri and Caron |355j . They considered the 
standard DMRG taking into account several low- lying target states (see Eq. ^ 
to construct the density matrix, weighted with the Boltzmann factor (/3 is the 
inverse temperature): 



With this method they performed reliable calculations of the magnetic sus- 
ceptibility of quantum spin chains with 5 = 1/2 and 3/2, showing excellent 
agreement with Bethe Ansatz exact results. They also calculated low temper- 
ature thermodynamical properties of the ID Kondo Lattice Model j356j and 
of organic conductors j357j . Zhang et al. j358j applied the same method in the 
study of a magnetic impurity embedded in a quantum spin chain. 

An alternative way of incorporating temperature, called the Transfer Ma- 
trix DMRG (TDMRG) stems from the Trotter-Suzuki expansion of the par- 
tition function of a one-dimensional model, turning it into a classical two- 
dimensional model, with the new axis corresponding to imaginary time (in- 
verse temperature). The extension of the DMRG method to classical sys- 
tems paved the way to study ID quantum systems at non zero tempera- 
ture ^linnnilSnilinniEnSKsee also a recent review in Ref. [2HS])- In this 
case the system is infinite and the finiteness is in the level of the Trotter ap- 
proximation. The free energy in the thermodynamic limit is determined by the 
largest eigenvalue of the transfer matrix. Being a high-temperature expansion, 
the method loses precision at low temperatures. 

Very nice results have been obtained for the dimerized, S = 1/2, XY model, 
where the specific heat was calculated involving an extremely small basis 
set |359) (m = 16), the agreement with the exact solution being much better 
in the case where the system has a substantial gap. It has also been used to 
calculate thermodynamic properties of the anisotropic S = 1/2 Heisenberg 
model, with relative errors for the spin susceptibility of less than 10~^ down 
to temperatures of the order of O.OIJ keeping m = 80 states |361j . A complete 
study of thermodynamic properties like magnetization, susceptibility, specific 
heat and temperature dependent correlation functions for the 5 = 1/2 and 
3/2 Heisenberg models was done in (363) . Other applications have been the 
calculation of the temperature dependence of the charge and spin gap in the 
Kondo insulator 364 , the calculation of thermodynamic properties of ferri- 
magnetic chains |,365j and spin ladders |136j , the study of impurity properties 
in spin chains fSBSirSBTj . frustrated quantum spin chains (368] . t-J j369j and 




(34) 
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spin ladders ! 37i)j and dimerized frustrated Heisenberg chains |371j . Recent 
studies include the thermodynamics and crossover phenomena in the ID t-J 
model using this method which led to very accurate results, as compared to 
the particular case of the supersymmetric limit where analytic Bethe-Ansatz 
results are available j372j and the calculation of the specific heat in doped 
anisotropic Hubbard ladders with charge-order instability j373j . 

Recently, new proposals for calculating finite-temperature properties arose 
based on the imaginary-time evolution of the wave function in the matrix- 
product ^22,.375..374, and ancillary states [22 ll376j formalisms to simulate evo- 
lutions of physical quantities in real and imaginary time and at finite tempera- 
ture. They are based on the so-defined projected entangled-pair states (PEPS) 
and the matrix product density operators (MPDO). The evolution in imagi- 
nary time to construct thermal and mixed states leads to a versatile DMRG 
algorithm that is restricted neither to large temperatures nor to homogeneous 
systems and opens the possiblity of simulating finite temperature, dissipation 
and decoherence effects. From the DMRG point of view, the PEPS or ancilla 
method is particularly appropriate. The Hilbert space is enlarged by introduc- 
ing auxiliary states (in the form of an extra chain coupled to the system for 
example) and then the imaginary-time evolution of a pure state in the enlarged 
space is calculated (this process is known as "purification" in quantum infor- 
mation theory HHI)- The ancilla sites serve as a thermodynamic bath which is 
traced out to obtain the thermodynamic behaviour of the original system. 



II 




Q Q Q Q 



System 



(b 6 6 



Environment 



Figure 15. Schematic representation of the ancilla approach used to calculate imaginary-time 
evolution. The real system corresponds to the lower sites (I) and the ancillary states are labelled by 
II. Full lines represent maximally entangled states. 



Let our original system be defined by states \si) and our auxiliary bath by 
Si) (Fig. I15() . Following the notation of |376|, we define the non-normalized 
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pure state 

|V(/?))=e-^^/2|V^(0)) (35) 

where the Hamiltonian H acts only on the real system (I). If /3 = (T ^ oo) 
this wave function is in a maximally entangled state between the original and 
auxiliary states: 



\m) = WY.\'^'^) (36) 



Si 



The partition function for arbitrary (3 is S = (^(/3)|'(/'(/3)) and the thermal 
average of an operator A is calculated as 

To obtain the temperature-dependent wave function, one can start with the 
maximally entangled state |V'(0)) and then evolve in imaginary time applying 
the Hamiltonian to the real system only. Within the DMRG framework, the 
blocks can be thought of as for example shown in Fig |15l where each pair 
\siSi) forms a supersite and the maximally entangled state is represented as a 
full line. The imaginary-time evolution is performed with any of the real-time 
algorithms (for example, the adaptive algorithm 295,290^ for nearest neigh- 
bour Hamiltonians). In |,S7Hj . the specific heat, magnetic susceptibility and the 
correlation length are calculated in spin chains with an accuracy comparable 
to the transfer- matrix DMRG (see Figilfi)). 



10.1 Finite-temperature dynamics 

In order to include temperature in the calculation of dynamical quantities, the 
Transfer Matrix RG described above (TMRG) |359li;-{m [i:S(i2] was extended to 
obtain imaginary time correlation functions |377p378p379) . After Fourier trans- 
formation in the imaginary time axis, analytic continuation from imaginary to 
real frequencies is done using maximum entropy (ME). The combination of the 
TMRG and ME is free from statistical errors and the negative sign problem of 
Monte Carlo methods but still has the extrapolation error of the analytic con- 
tinuation. Since we are dealing with the transfer matrix, the thermodynamic 
limit (infinite system size) can be discussed directly without extrapolations. 
However, in the present scheme, only local quantities can be calculated. 

A systematic investigation of local spectral functions is done in Ref. j379j 
for the anisotropic Heisenberg antiferromagnetic chain. The authors obtain 



New Trends in Density Matrix Renormalization 



41 




Figure 16. Specific heat and magnetic susceptibility of the Haldane chain (S = 1) with 64 sites 
calculated using the Trotter-Suzuki time-evolution algorithm. Also shown are the results using 
Transfer Matrix DMRG (reprinted from Ref. 18761 with permission). 

good qualitative results especially for high temperatures but a quantitative 
description of peaks and gaps are beyond the method, due to the severe in- 
trinsic limitation of the analytic continuation. This method was also applied 
with great success to the ID Kondo insulator |378j where the temperature 
dependence of the local density of states and local dynamic spin and charge 
correlation functions were calculated. A modfication of this method to avoid 
the use of the ill-posed analytic continuation was done in |38n) by considering 
a path-integral approach to calculate real-time correlation functions at finite 
temperatures. 



11 Dynamical Mean Field Theory using DMRG 

The Dynamical Mean Field Theory (DMFT) has become one of the basic 
methods to calculate realistic electronic band structure in strongly correlated 
systems |381) . At the heart of the DMFT method is the solution of an associ- 
ated quantum impurity model where the environment of the impurity has to 
be determined self-consistently. Therefore the ability to obtain reliable DMFT 
solutions of lattice model Hamiltonians relies directly on the ability to solve 
quantum impurity models. Among the a priori exact numerical algorithms 
available we find the Hirsch-Fye Quantum Monte Carlo '382 , 383| method and 
Wilson's Numerical Renormalization Group (NRG) (8k384..385) . While the for- 
mer, a finite-temperature method, is very stable and accurate at the Matsubara 
frequencies, its main drawback is the access to real frequency quantities for the 
calculation of spectral functions which requires less controlled techniques for 
the analytic continuation of the Green functions. The second method can be 
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formulated both at T=0 and finite (small) T and provides extremely accurate 
results at very small frequencies, at the expense of a less accurate description 
of the high energy features. 

In order to overcome the difficulties encountered by these other methods, we 
will show that DMRG can be used very reliably to solve the related impurity 
problem within DMFT |386j . By using the DMRG to solve the impurity, no 
a priori approximations are made and the method provides equally reliable 
solutions for both gapless and gapful phases. More significantly, it provides 
accurate estimates for the distributions of spectral intensities of high frequency 
features such as the Hubbard bands, that are of main relevance for analysis of 
X-ray photoemission and optical conductivity experiments. 

We will now very briefly describe the method applied to the Mott transition 
in the Hubbard model. The Hamiltonian of the Hubbard model is defined by 
Eq. 1151 Applying DMFT to this model leads to a mapping of the original lat- 
tice model onto an associated quantum impurity problem in a self-consistent 
bath. In the particular case of the Hubbard model, the associated impurity 
problem is the single impurity Anderson model (SIAM), where the hybridiza- 
tion function A(a;), which in the usual SIAM is a flat density of states of 
the conduction electrons, is now to be determined self-consistently. More pre- 
cisely, for the Hamiltonian ()15() defined on a Bethe lattice of coordination d, 
one takes the limit of large d and exactly maps the model onto a SIAM im- 
purity problem with the requirement that A(a;) = t^G{ijj), where G(a;) is the 
impurity Green's function. At the self-consistent point G{lli) coincides with 
the local Green's function of the original lattice model |387j . We take the 
half-bandwidth of the non interacting model as unit of energy, t = 1/2. 

The Green's function of the impurity problem is an important quantity in 
this algorithm: Go(a;) = l/(w + /x - A(u;)) = l/{uj + n - t^G{uj)). Thus, 
to implement the new algorithm we shall consider j3S8|l,S89] a general repre- 
sentation of the hybridization function in terms of continued fractions that 
define a parametrization of A(lj) in terms of a set of real and positive coef- 
ficients. Since it is essentially a Green's function, A{z) can be decomposed 
into "particle" and "hole" contributions as A(z) = A''{z) + A.^{z) with 
A>{z) = t^{gs\cj^^^::^c^gs) and A<{z) = t^{gs\c^ j^^ff^c\gs) for a 
given Hamiltonian, H with ground-state energy Eq. By standard Lanczos tech- 
nique, H can be tridiagonalized and the functions A'^{z) and A^(z) can be 
expressed in terms of respective continued fractions jHS] • As first implemented 
in Ref. |388[I389] . each continued fraction can be represented by a chain of 
auxiliary atomic sites whose energies and hopping amplitudes are given by the 
continued fraction diagonal and off-diagonal coefficients respectively. 

As a result of the self-consistency condition, the two chains representing the 
hybridization, are "attached" to the right and left of an atomic site to obtain a 
new SIAM Hamiltonian, H. In fact Go (-2;) constitutes the local Green's function 
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of the site plus chain system. The algorithm in Ref. [HHH'SS^, basically consists 
in switching on the local Coulomb interaction at the impurity site of the SIAM 
Hamiltonian and using the Lanczos technique to re-obtain A{z), iterating the 
procedure until the set of continued fractions coefficients converges. By using 
the DMRG, the number of auxiliary sites that can be used in the hybridization 
chains is much larger than in the exact diagonalization scheme, leading to 
more accurate results |39U| . An alternative way of using DMRG to solve the 
impurity problem, which does not rely on the continued fraction expansion 
was developed in |391j . 

The SIAM Hamiltonian therefore reads 

Nc Nc-1 
a,a=-Nc;ajtO o-,a=-(7Vc-l);a7^0,-l 

+ Yl bo{cUaa + h.c.)+U{n^ -h{ni-^) (38) 



with Co- being the destruction operator at the impurity site, and Caa being the 
destruction operator at the a site of the hybridization chain of 2Nc sites. The 
set of parameters {aa,ba} are directly obtained from the coefficients of the 
continued fraction representations of A{z) by the procedure just described. 

In Fig llTI we show the DMFT+DMRG results (solid lines) for the density 
of states (DOS) for several values of increasing interaction U. The results are 
compared to the Iterated Perturbation Theory (IPT) results (dashed lines) 
P2,393 . 

Results for the imaginary part of the Green's funtions on the Matsubara 
frequencies match the precise Monte Carlo solutions at low temperatures. We 
also obtain accurate values for the two distinct critical values of the interaction 
Uci == 2.39 lb 0.02 and Uc2 = 3.0 it 0.2. To this end, reaching larger system 
sizes turned out to be important in order to perform proper extrapolations 
and overcome finite-size effects. 

As a conclusion, using DMRG as the impurity-solver of the DMFT, increases 
its performance. Large systems can be considered and accurate values of the 
critical interactions are obtained in agreement with NRG predictions allowing 
for a non-trivial test of the accuracy of this method. In contrast with NRG, 
however, this new algorithm deals with all energy scales on equal footing which 
allowed us to find interesting substructure in the Hubbard bands of the cor- 
related metallic state. The ability of the new algorithm to directly deal with 
the high energy scales is a very important feature which is relevant for the 
interpretation of high resolution photoemission spectroscopies. |381| In addi- 
tion, with this method, realistic band-structure calculations of systems with a 
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Figure 17. Density of states of the half-filled Hubbard Model. We also show the IPT results 

(dashed lines) (see text). 



larger number of degrees of freedom and other dopings away from half-filling 
can be handled ??. 



12 Summary and outlook 

The DMRG method has proven to be a very reliable and versatile numerical 

method that can be applied to a broad spectrum of problems in physics in a 
variety of fields such as Condensed Matter, Statistical Mechanics, High Energy 
and Nuclear Physics, Quantum Chemistry and Quantum Information Theory. 
It is nowadays recognized as one of the most accurate and efficient numerical 
techniques. We have reviewed the most recent developments to DMRG among 
which are the very accurate real-time applications, the use of quantum infor- 
mation concepts which give an interesting perspective and the possibility of 
using it as the impurity-solver in DMFT calculations. These new improvements 
have triggered a great activity and numerous papers have been published ap- 
plying these new techniques. 

However, the DMRG is still evolving and presents very interesting potential 
new applications. Among the fields where DMRG has still to be exploited are 
the time-dependent analysis of quantum and dissipative systems, the appli- 
cation in DMFT band-structure calculations of more complex and realistic 
systems and quantum chemistry. The connection of DMRG with quantum in- 
formation theory paves the way to an important new field of application and, 
very possibly, to dimensions higher than one. 
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In view of the fact that numerical simulations constitute an essential tool in 
physics and quantum chemistry, a group of experts have started an interesting 
initiative called the ALPS (Algorithms and Libraries for Physics Simulations) 
project |395j . It consists of an open source of optimized software aimed at 
studying strongly correlated systems such as quantum magnets, strongly cor- 
related fermionic systems and lattice bosons where methods like classical and 
quantum Monte Carlo, exact diagonalization and DMRG are included. 
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